home *** CD-ROM | disk | FTP | other *** search
/ APDL Eductation Resources / APDL Eductation Resources.iso / programs / electronic / rlab / TestMatrix / clement_r < prev    next >
Encoding:
Text File  |  1994-12-20  |  1.5 KB  |  55 lines
  1. //-------------------------------------------------------------------//
  2.  
  3. // Synopsis:    Clement matrix - tridiagonal with zero diagonal entries.
  4.  
  5. // Syntax:    A = clement ( N , K )
  6.  
  7. // Description:
  8.  
  9. //    A is a tridiagonal matrix with zero diagonal entries and known
  10. //    eigenvalues.  It is singular if N is odd.  About 64 percent of
  11. //    the entries of the inverse are zero.  The eigenvalues  are
  12. //    plus and minus the numbers N-1, N-3, N-5, ..., (1 or 0).
  13.  
  14. //    For K = 0 (the default) the matrix is unsymmetric, while for 
  15. //    K = 1 it is symmetric. 
  16.  
  17. //    clement(N, 1) is diagonally similar to clement(N).
  18.  
  19. //    Similar properties hold for tridiag(X,Y,Z) where 
  20. //    Y = zeros(N,1). The eigenvalues still come in plus/minus pairs
  21. //    but they are not known explicitly.
  22.  
  23. //      References:
  24. //        P.A. Clement, A class of triple-diagonal matrices for test
  25. //             purposes, SIAM Review, 1 (1959), pp. 50-52.
  26. //        O. Taussky and J. Todd, Another look at a matrix of Mark Kac,
  27. //             Linear Algebra and Appl., 150 (1991), pp. 341-360.
  28.  
  29. //    This file is a translation of clement.m from version 2.0 of
  30. //    "The Test Matrix Toolbox for Matlab", described in Numerical
  31. //    Analysis Report No. 237, December 1993, by N. J. Higham.
  32.  
  33. //-------------------------------------------------------------------//
  34.  
  35. clement = function ( n , k )
  36. {
  37.   local (n, k)
  38.  
  39.   if (!exist (k)) { k = 0; }
  40.  
  41.   n = n-1;
  42.   x = n:1:-1;
  43.   z = 1:n;
  44.  
  45.   if (k == 0)
  46.   {
  47.     A = diag(x, -1) + diag(z, 1);
  48.   else
  49.     y = sqrt(x.*z);
  50.     A = diag(y, -1) + diag(y, 1);
  51.   }
  52.  
  53.   return A;
  54. };
  55.